# load libraries
library(tidyverse)
library(psych)
library(langcog) # source: https://github.com/langcog/langcog
library(RColorBrewer)
library(plotly)
library(lubridate)
library(rms)
library(lme4)
library(brms)
library(parallel)
# clear workspace
# rm(list = ls(all = T))
graphics.off()
# set number of cores for brms
n_cores <- detectCores()
# make rounding function
round2 <- function(x) {format(round(x, 2), nsmall = 2)}
# make cleanup function
cleanup <- function(datasource, age_group) {
if(grepl("adult", age_group)) {
# set target dataset
if(datasource == "study 1"){d <- d_raw_study1}
if(datasource == "study 1b"){d <- d_raw_study1b}
if(datasource == "study 1c"){d <- d_raw_study1c}
# enact exclusionary criteria
d_clean_1 <- d
# recode background and demographic variables
d_clean <- d_clean_1 %>%
mutate( # deal with study number
study = factor(study)) %>%
mutate( # deal with race
race_cat2 = factor(sub(" +$", "", ethnicity)),
race_cat3 = factor(ifelse(grepl(" ", race_cat2) == T, "multiracial",
as.character(race_cat2)))) %>%
dplyr::select(study, subid:country_selfrep, age_group, race_cat3) %>%
rename(race_cat = race_cat3) %>%
mutate( # deal with religion (note: only dealing with childhood religion for now)
religion_cat2 = factor(sub(" +$", "", religionChild)),
religion_cat3 = factor(ifelse(grepl(" ", religion_cat2) == T,
"multireligious",
as.character(religion_cat2)))) %>%
dplyr::select(study:race_cat, religion_cat3) %>%
rename(religion_cat = religion_cat3)
# remove extraneous dfs and variables
rm(d, d_clean_1)
}
if(grepl("child", age_group)) {
# set target dataset
if(datasource == "study 2"){d <- d_raw_study2}
# if(datasource == "study 3"){d <- d_raw_study3}
# if(datasource == "study 4"){d <- d_raw_study4}
# recode background and demographic variables
d_clean_2 <- d %>%
mutate( # deal with study number
study = factor(study),
responseNum = ifelse(!is.na(responseNum), responseNum,
ifelse(response == "no", 0,
ifelse(response == "kinda", 0.5,
ifelse(response == "yes", 1, NA)))))
# NOTE: need to reconcile race/ethnicity at some point...
# NOTE: need to deal with gender at some point...
d_clean <- d_clean_2
# remove extraneous dfs and variables
rm(d, d_clean_2)
}
# remove outliers if desired
if(chosenOutlierHandling == "remove") {
d_clean <- d_clean %>%
gather(capacity, score, happy:pride) %>%
group_by(character, capacity) %>%
filter(!score %in% boxplot.stats(score, 2.5)$out) %>%
spread(capacity, score) %>%
arrange(character, subid)
}
# filter characters if desired
if(is.element("none", chosenExclude)) {} else {
d_clean <- d_clean %>%
filter(!character %in% chosenExclude)
}
# filter items if desired
if(is.element("none", chosenExcludeItem)) {} else {
d_clean <- d_clean %>%
dplyr::filter(!capacity %in% chosenExcludeItem)
}
# drop trials <250 ms
d_clean <- d_clean %>%
filter(rt >= 250 | is.na(rt))
# center response variable
if(datasource == "study 1b") {
d_clean <- d_clean %>%
mutate(responseNumC = responseNum - 4)
} else {
d_clean <- d_clean %>%
mutate(responseNumC = responseNum - 0.5)
}
# rename character name variables
if("charName" %in% names(d_clean)) {
d_clean <- d_clean %>% rename(character = charName)
}
# cleanup
d_clean <- d_clean %>%
filter(!is.na(subid), !is.na(character), !is.na(capacity))
# return cleaned dataset
return(d_clean)
}
# make function for stripping dataframes for dimension reducation
makeDRDF <- function(datasource, chosenCondition) {
# set target dataset
if(datasource == "study 1"){d <- d1}
if(datasource == "study 1b"){d <- d1b}
if(datasource == "study 1c"){d <- d1c}
if(datasource == "study 2"){d <- d2}
# filter by character if specified
if(chosenCondition %in% c("beetle", "robot")) {
d <- d %>% filter(character == chosenCondition)
}
# make stripped dataframe for dimension reducation analyses
d_strip <- d %>%
filter(!is.na(character), !is.na(subid), !is.na(capacity), capacity != "") %>%
mutate(subid = paste(character, subid, sep = "_")) %>%
select(subid, capacity, responseNum) %>%
spread(capacity, responseNum) %>%
remove_rownames() %>%
column_to_rownames(var = "subid")
# return stripped dataframe
return(d_strip)
}
# make demographics functions
demoSampleSize <- function(datasource) {
# set target dataset
if(datasource == "study 1"){d <- d1}
if(datasource == "study 1b"){d <- d1b}
if(datasource == "study 1c"){d <- d1c}
if(datasource == "study 2"){d <- d2}
# get distinct subids
sample_size <- d %>% distinct(subid, character) %>% count(character) %>% data.frame()
# add total sample size
sample_size <- rbind(sample_size %>% mutate(character = as.character(character)),
c(character = "all", n = d %>% distinct(subid) %>% count() %>% as.numeric()))
# return dataframe
return(sample_size)
}
demoDuration <- function(datasource) {
# set target dataset
if(datasource == "study 1"){d <- d1}
if(datasource == "study 1b"){d <- d1b}
if(datasource == "study 1c"){d <- d1c}
if(datasource == "study 2"){d <- d2 %>% rename(duration = sessionDuration)}
# get sample size per character
duration <- d %>%
distinct(subid, character, duration) %>%
mutate(duration = as.numeric(duration)) %>%
group_by(character) %>%
summarise(min_duration = min(duration, na.rm = T),
max_duration = max(duration, na.rm = T),
median_duration = median(duration, na.rm = T),
mean_duration = mean(duration, na.rm = T),
sd_duration = sd(duration, na.rm = T))
# add total duration
all <- d %>%
distinct(subid, character, duration) %>%
mutate(duration = as.numeric(duration)) %>%
summarise(min_duration = min(duration, na.rm = T),
max_duration = max(duration, na.rm = T),
median_duration = median(duration, na.rm = T),
mean_duration = mean(duration, na.rm = T),
sd_duration = sd(duration, na.rm = T)) %>%
mutate(character = "all")
duration <- rbind(duration, all) # not sure why full_join doesn't work
# return dataframe
return(duration)
}
demoAge <- function(datasource) {
# set target dataset
if(datasource == "study 1"){d <- d1}
if(datasource == "study 1b"){d <- d1b}
if(datasource == "study 1c"){d <- d1c}
if(datasource == "study 2"){d <- d2}
# get sample size per character
age <- d %>%
distinct(subid, character, age) %>%
mutate(age = as.numeric(age)) %>%
group_by(character) %>%
summarise(min_age = min(age, na.rm = T),
max_age = max(age, na.rm = T),
median_age = median(age, na.rm = T),
mean_age = mean(age, na.rm = T),
sd_age = sd(age, na.rm = T))
# add total age
all <- d %>%
distinct(subid, character, age) %>%
mutate(age = as.numeric(age)) %>%
summarise(min_age = min(age, na.rm = T),
max_age = max(age, na.rm = T),
median_age = median(age, na.rm = T),
mean_age = mean(age, na.rm = T),
sd_age = sd(age, na.rm = T)) %>%
mutate(character = "all")
age <- full_join(age, all)
# return dataframe
return(age)
}
demoGender <- function(datasource) {
# set target dataset
if(datasource == "study 1"){d <- d1}
if(datasource == "study 1b"){d <- d1b}
if(datasource == "study 1c"){d <- d1c}
if(datasource == "study 2"){d <- d2}
# get gender per character and overall
gender <- data.frame(addmargins(with(d %>% distinct(subid, character, gender),
table(character, gender)))) %>%
filter(gender != "Sum") %>%
rename(n = Freq)
gender <- gender %>%
mutate(character = factor(ifelse(character == "Sum",
"all", as.character(character)),
levels = c("beetle", "robot", "all"))) %>%
arrange(character, gender) %>%
spread(gender, n)
# return dataframe
return(gender)
}
demoRace <- function(datasource) {
# set target dataset
if(datasource == "study 1"){d <- d1}
if(datasource == "study 1b"){d <- d1b}
if(datasource == "study 1c"){d <- d1c}
if(datasource == "study 2"){d <- d2 %>%
mutate(ethnicity = as.character(ethnicity),
white = grepl("white", ethnicity) | grepl("cauc", ethnicity) |
grepl("euro", ethnicity) | grepl("australia", ethnicity),
black = grepl("black", ethnicity) | grepl("africa", ethnicity),
latinx = grepl("hisp", ethnicity) | grepl("latin", ethnicity),
E_asian = (grepl("east asian", ethnicity) & grepl("southeast", ethnicity) == F) |
grepl("chin", ethnicity) | grepl("korea", ethnicity)
| grepl("\\(asian\\)", ethnicity) |
grepl("\\/asian", ethnicity),
SSE_asian = grepl("south or southeast", ethnicity) |
(grepl("india", ethnicity) & grepl("native", ethnicity) == F) |
grepl("vietnam", ethnicity) | grepl("pakistan", ethnicity),
nat_am = grepl("native american", ethnicity),
middle_east = grepl("middle", ethnicity)) %>%
mutate(race_cat = ifelse(white + black + latinx + E_asian + SSE_asian +
nat_am + middle_east > 1,
"multiracial",
ifelse(white, "white",
ifelse(black, "black",
ifelse(latinx, "latinx",
ifelse(E_asian, "east_asian",
ifelse(SSE_asian, "south_southeast_asian",
ifelse(nat_am, "native_american",
ifelse(middle_east, "middle_eastern",
NA)))))))))
}
# get race per character and overall
race <- data.frame(addmargins(with(d %>% distinct(subid, character, race_cat),
table(character, race_cat)))) %>%
filter(race_cat != "Sum") %>%
rename(n = Freq)
race <- race %>%
mutate(character = factor(ifelse(character == "Sum",
"all", as.character(character)))) %>%
arrange(character, race_cat) %>%
spread(race_cat, n)
# return dataframe
return(race)
}
# plotting functions
makeFacetLabs <- function(df_plotting) {
facet_labels <- array()
df_plotting <- df_plotting %>% mutate(character = factor(character))
for(i in 1:length(levels(df_plotting$character))) {
df <- df_plotting %>% filter(character == levels(df_plotting$character)[i]) %>%
select(character, n) %>% unique()
facet_labels[i] <- paste0(df$character, " (n = ", df$n, ")")
}
names(facet_labels) <- levels(df_plotting$character)
return(facet_labels)
}
# remove outliers?
chosenOutlierHandling <- "keep"
# chosenOutlierHandling <- "remove"
# exclude any conditions (characters)?
chosenExclude <- "none"
# chosenExclude <- c("stapler", "car", "computer")
# exclude any items (mental capacities)?
# chosenExcludeItem <- "none"
# chosenExcludeItem <- "computations"
chosenExcludeItem <- c("metal", "on_off")
# NOTE: always choose minimal residual (fm = "minres") instead of ML because of non-normality
# for EFAs, what kind of correlation?
chosenCorType <- "cor" # pearson correlation
# chosenCorType <- "poly" # polychoric correlation
# for EFAs, what kind of rotation?
chosenRotType <- "oblimin" # oblimin rotation
# chosenRotType <- "oblimin" # oblimin rotation
# chosenRotType <- "none" # no rotation
data.frame("conditionsExcluded" = chosenExclude,
"outlierHandling" = chosenOutlierHandling,
"EFA_correlation" = chosenCorType,
"EFA_rotation" = chosenRotType)
# study 1 (2016-07-06, adults, 2 conditions, 3-point scale, "decide what to do" and "make plans")
d_raw_study1 <- read.csv("/Users/kweisman/Documents/Research (Stanford)/Projects/Dimkid/dimkid/data/adults/us_run-01_2016-06-05_anonymized.csv") %>%
mutate(study = "study 1", age_group = "adults") %>% select(-X)
# study 1b (2017-07-19, adults, 2 conditions, 7-point scale, "decide what to do" and "make plans")
d_raw_study1b <- read.csv("/Users/kweisman/Documents/Research (Stanford)/Projects/Dimkid/dimkid/data/adults/us_run-02_2016-07-19_anonymized.csv") %>%
mutate(study = "study 1b", age_group = "adults") %>% select(-X)
# study 1c (2016-12-08, adults, 2 conditions, 3-point scale, "have free will" and "have intentions")
d_raw_study1c <- read.csv("/Users/kweisman/Documents/Research (Stanford)/Projects/Dimkid/dimkid/data/adults/us_run-03_2016-12-08_anonymized.csv") %>%
mutate(study = "study 1c", age_group = "adults") %>% select(-X)
# study 2 (June - December 2016, 7-9yo, 2 conditions, 3-point-scale, "decide what to do" and "make plans")
d_raw_study2 <- read.csv("/Users/kweisman/Documents/Research (Stanford)/Projects/Dimkid/dimkid/data/children/run-01_2017-07-24_anonymized.csv") %>%
mutate(study = "study 2", age_group = "children_79") %>% select(-X)
# clean up datasets
d1 <- cleanup("study 1", "adults")
d1b <- cleanup("study 1", "adults")
d1c <- cleanup("study 1", "adults")
d2 <- cleanup("study 2", "children")
# tweak by hand
d2 <- d2 %>%
filter(!is.na(age)) %>%
filter(age >= 7, age < 10) %>%
filter(character != "elephant")
# fast rt
d_raw_study1 %>%
full_join(d_raw_study2) %>%
filter(!is.na(rt)) %>%
mutate(rt_cat = ifelse(rt < 250, "<250ms", ">=250ms")) %>%
count(age_group, rt_cat) %>%
group_by(age_group) %>%
mutate(prop = round(n/sum(n), 2))
Joining, by = c("run", "subid", "age", "gender", "ethnicity", "trialNum", "bgColor", "capacity", "capWording", "hoverTime", "rt", "response", "responseNum", "study", "age_group")
Column `run` joining factors with different levels, coercing to character vectorColumn `subid` joining factors with different levels, coercing to character vectorColumn `gender` joining factors with different levels, coercing to character vectorColumn `ethnicity` joining factors with different levels, coercing to character vectorColumn `capWording` joining factors with different levels, coercing to character vectorColumn `response` joining factors with different levels, coercing to character vector
# skipped trials
d_raw_study1 %>%
full_join(d_raw_study2) %>%
filter((age >= 7 & age < 10) | (age >= 18)) %>%
filter(!is.na(response), !response %in% c("bail", "skip")) %>%
count(age_group, subid) %>%
filter(n != 40) %>%
mutate(skipped = 40 - n) %>%
group_by(age_group) %>%
summarise(total_skipped = sum(skipped))
Joining, by = c("run", "subid", "age", "gender", "ethnicity", "trialNum", "bgColor", "capacity", "capWording", "hoverTime", "rt", "response", "responseNum", "study", "age_group")
Column `run` joining factors with different levels, coercing to character vectorColumn `subid` joining factors with different levels, coercing to character vectorColumn `gender` joining factors with different levels, coercing to character vectorColumn `ethnicity` joining factors with different levels, coercing to character vectorColumn `capWording` joining factors with different levels, coercing to character vectorColumn `response` joining factors with different levels, coercing to character vector
# total mising
d1 %>%
full_join(d2) %>%
filter((age >= 7 & age < 10) | (age >= 18)) %>%
filter(!is.na(response), !response %in% c("bail", "skip")) %>%
count(age_group, subid) %>%
filter(n != 40) %>%
mutate(missing = 40 - n) %>%
group_by(age_group) %>%
summarise(total_missing = sum(missing),
prop_missing = round(total_missing/(40*200), 3))
Joining, by = c("study", "subid", "character", "age", "gender", "ethnicity", "trialNum", "bgColor", "capacity", "capWording", "hoverTime", "rt", "response", "responseNum", "age_group", "responseNumC")
Column `study` joining factors with different levels, coercing to character vectorColumn `subid` joining factors with different levels, coercing to character vectorColumn `character` joining factors with different levels, coercing to character vectorColumn `gender` joining factors with different levels, coercing to character vectorColumn `ethnicity` joining factors with different levels, coercing to character vectorColumn `capWording` joining factors with different levels, coercing to character vectorColumn `response` joining factors with different levels, coercing to character vector
# make dataframes for s1
# d1_beetle <- makeDRDF("study 1", "beetle")
# d1_robot <- makeDRDF("study 1", "robot")
d1_all <- makeDRDF("study 1", "all")
# make dataframes for follow-up studies to s1
d1b_all <- makeDRDF("study 1b", "all")
d1c_all <- makeDRDF("study 1c", "all")
# make dataframes for study 2
# d2_beetle <- makeDRDF("study 2", "beetle")
# d2_robot <- makeDRDF("study 2", "robot")
d2_all <- makeDRDF("study 2", "all")
For all studies we conduct exploratory factor analyses using Pearson correlations to find minimum residual solutions.
For each study, we first examine maximal unrotated and rotated solutions. To determine the maximum number of factors to extract, we use the following rule of thumb: With \(p\) observations per participant, we can extract a maximum of \(k\) factors, where \((p-k)*2 > p+k\), i.e., \(k < p/3\). Thus, with 40 mental capacity items, we can extract a maximum of 13 factors.
To determine how many factors to retain, we use the following preset retention criteria, considering the unrotated maximal solution (unless otherwise noted):
We then examine and interpret oblimin-rotated solutions, extracting only the number of factors that meet these criteria.
Study information:
# make demographics tables
demoSampleSize("study 1")
demoDuration("study 1")
demoAge("study 1")
Joining, by = c("character", "min_age", "max_age", "median_age", "mean_age", "sd_age")
Column `character` joining factor and character vector, coercing into character vector
demoGender("study 1")
demoRace("study 1") %>%
filter(character == "all") %>%
gather(ethnicity, count, -character) %>%
mutate(prop = round(count/sum(count), 2)) %>%
arrange(desc(prop))
# gender
with(d1 %>%
distinct(subid, character, gender) %>%
filter(!is.na(gender)) %>%
mutate(character = factor(character)),
table(character, gender)) %>%
summary()
Number of cases in table: 200
Number of factors: 2
Test for independence of all factors:
Chisq = 1.9863, df = 2, p-value = 0.3704
Chi-squared approximation may be incorrect
# age
with(d1 %>%
distinct(subid, character, age),
t.test(age ~ character))
Welch Two Sample t-test
data: age by character
t = 0.26048, df = 192.06, p-value = 0.7948
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-2.625759 3.424833
sample estimates:
mean in group beetle mean in group robot
33.15464 32.75510
# examine scree plot
# fa.parallel(d1_all)
# run EFA without rotation with N factors
efa_d1_all_unrotated <- fa(d1_all, 13, rotate = "none",
cor = chosenCorType, fm = "minres")
fa.sort(efa_d1_all_unrotated)
Factor Analysis using method = minres
Call: fa(r = d1_all, nfactors = 13, rotate = "none", fm = "minres",
cor = chosenCorType)
Standardized loadings (pattern matrix) based upon correlation matrix
MR1 MR2 MR3 MR4 MR5 MR6 MR7 MR8 MR9 MR10 MR11
happy 0.76 0.00 -0.33 -0.08 -0.22 0.08 -0.15 -0.03 -0.08 -0.04 -0.08
joy 0.76 0.01 -0.39 0.10 -0.17 -0.01 -0.08 -0.07 0.05 0.00 -0.03
love 0.75 0.11 -0.28 0.09 -0.07 -0.09 -0.07 -0.14 -0.03 -0.08 0.00
depressed 0.74 0.04 -0.37 0.04 -0.14 0.04 -0.17 -0.04 -0.21 -0.05 -0.02
fear 0.72 -0.39 0.14 0.03 -0.18 0.07 0.09 0.12 -0.03 0.13 -0.04
safe 0.71 -0.29 0.21 -0.12 -0.03 -0.06 0.04 0.06 -0.03 0.09 -0.13
tired 0.69 -0.34 0.23 0.07 -0.06 0.05 0.06 0.08 0.07 -0.07 0.09
pleasure 0.69 -0.23 -0.07 0.15 -0.15 -0.08 -0.15 0.06 0.16 -0.06 0.21
calm 0.68 -0.17 0.01 -0.08 -0.08 0.04 0.08 0.10 -0.04 0.02 -0.01
pride 0.68 0.18 -0.42 0.08 0.04 0.04 0.01 -0.08 -0.14 0.05 -0.06
desires 0.66 -0.17 0.10 -0.02 -0.03 -0.02 0.35 -0.48 0.13 -0.14 -0.08
angry 0.65 -0.04 -0.11 -0.03 0.05 0.22 -0.02 0.04 -0.10 0.17 0.02
nauseated 0.65 -0.32 0.14 0.08 -0.16 0.05 0.01 -0.06 0.14 0.03 0.03
pain 0.63 -0.52 0.19 -0.01 -0.04 0.11 0.11 0.19 0.05 -0.01 0.09
disrespected 0.63 0.06 -0.35 0.16 0.07 0.07 -0.13 -0.03 -0.07 -0.03 -0.07
guilt 0.62 0.14 -0.41 0.21 0.43 0.14 0.04 0.02 0.08 0.00 -0.03
thoughts 0.55 0.18 0.10 -0.37 -0.01 0.04 -0.04 0.10 -0.09 -0.15 0.19
embarrassed 0.52 0.14 -0.40 0.19 0.48 0.18 0.11 0.11 0.12 -0.03 0.09
odors 0.49 -0.35 0.37 0.05 0.15 -0.09 -0.03 0.10 0.07 0.01 -0.03
beliefs 0.48 0.40 -0.16 -0.14 0.04 -0.38 0.18 0.09 0.07 0.04 -0.06
self_aware 0.46 0.18 0.22 -0.30 0.09 0.00 0.12 0.13 -0.20 -0.28 -0.21
personality 0.44 0.36 -0.19 -0.13 -0.03 -0.27 0.01 0.00 0.24 0.20 -0.06
self_restraint 0.43 0.35 -0.05 -0.15 0.04 -0.19 -0.08 0.07 0.00 0.02 0.17
goal 0.41 0.21 0.19 -0.11 0.07 -0.11 0.18 -0.08 0.02 -0.03 0.07
choices 0.37 0.34 0.36 -0.20 0.06 0.09 -0.25 -0.20 0.12 -0.05 0.17
computations -0.33 0.82 -0.07 0.14 -0.03 0.00 -0.02 0.04 0.03 0.00 0.03
recognizing 0.10 0.76 0.12 0.13 -0.21 0.15 0.11 0.11 0.15 -0.02 0.00
hungry 0.55 -0.71 0.22 -0.03 0.08 -0.02 0.06 0.02 0.03 0.03 0.05
remembering 0.14 0.66 0.16 0.10 -0.15 0.13 0.01 -0.05 0.03 -0.15 0.04
communicating 0.11 0.62 0.18 0.14 -0.18 0.30 0.11 0.18 0.15 -0.05 -0.08
intentions 0.19 0.62 0.02 -0.16 0.01 0.07 0.26 0.00 -0.20 0.00 0.09
morality 0.31 0.50 -0.13 0.02 -0.07 -0.19 0.08 0.05 -0.04 0.00 0.18
reasoning 0.34 0.44 0.31 -0.16 0.01 0.21 -0.11 0.01 0.08 0.20 -0.21
emo_recog 0.37 0.39 -0.10 -0.10 0.01 -0.27 -0.09 0.10 0.04 0.08 -0.14
seeing 0.33 0.15 0.50 0.28 0.08 -0.03 -0.07 -0.01 -0.11 -0.01 0.07
depth 0.26 0.27 0.48 0.28 0.12 -0.09 -0.16 -0.09 0.04 -0.08 -0.10
temperature 0.30 0.19 0.46 0.40 0.05 -0.22 0.04 0.00 -0.26 0.12 0.09
conscious 0.44 0.10 0.44 -0.11 0.17 -0.11 -0.24 0.06 0.08 -0.19 -0.15
sounds 0.27 0.20 0.42 0.38 -0.06 -0.02 0.06 -0.07 -0.11 0.12 -0.07
free_will 0.31 0.30 0.32 -0.40 0.15 0.22 -0.07 -0.19 -0.05 0.27 0.10
MR12 MR13 h2 u2 com
happy -0.04 -0.08 0.80 0.20 1.8
joy 0.09 -0.01 0.79 0.21 1.8
love 0.08 0.07 0.71 0.29 1.6
depressed 0.06 -0.03 0.79 0.21 2.0
fear -0.16 0.15 0.82 0.18 2.3
safe -0.04 0.07 0.68 0.32 1.8
tired -0.11 -0.05 0.70 0.30 2.0
pleasure 0.09 0.18 0.72 0.28 2.3
calm 0.05 -0.07 0.54 0.46 1.3
pride -0.16 0.01 0.74 0.26 2.2
desires 0.02 -0.07 0.87 0.13 2.9
angry -0.20 -0.09 0.58 0.42 1.8
nauseated 0.04 0.05 0.61 0.39 2.0
pain -0.01 -0.14 0.80 0.20 2.7
disrespected -0.02 0.02 0.59 0.41 2.0
guilt -0.03 -0.05 0.84 0.16 3.3
thoughts 0.09 -0.05 0.57 0.43 2.9
embarrassed 0.12 0.07 0.82 0.18 4.3
odors 0.11 0.05 0.56 0.44 3.4
beliefs -0.09 -0.04 0.64 0.36 4.0
self_aware 0.15 -0.07 0.60 0.40 5.6
personality 0.09 -0.17 0.59 0.41 5.2
self_restraint -0.04 -0.14 0.43 0.57 3.5
goal -0.14 0.15 0.36 0.64 4.0
choices -0.07 -0.05 0.59 0.41 6.1
computations -0.06 -0.03 0.81 0.19 1.4
recognizing -0.02 -0.12 0.75 0.25 1.7
hungry 0.07 0.01 0.88 0.12 2.2
remembering 0.09 0.10 0.58 0.42 1.7
communicating -0.01 -0.03 0.65 0.35 2.6
intentions -0.01 0.18 0.60 0.40 2.3
morality -0.10 0.10 0.46 0.54 2.9
reasoning 0.09 0.12 0.60 0.40 5.4
emo_recog 0.14 0.11 0.46 0.54 4.3
seeing 0.09 -0.08 0.50 0.50 3.1
depth -0.13 -0.02 0.53 0.47 4.1
temperature 0.06 -0.09 0.65 0.35 4.8
conscious -0.21 0.05 0.63 0.37 4.6
sounds 0.08 0.01 0.48 0.52 4.0
free_will 0.08 0.00 0.65 0.35 6.8
MR1 MR2 MR3 MR4 MR5 MR6 MR7 MR8 MR9 MR10 MR11 MR12
SS loadings 11.03 5.51 3.15 1.29 0.87 0.84 0.64 0.56 0.51 0.46 0.43 0.38
Proportion Var 0.28 0.14 0.08 0.03 0.02 0.02 0.02 0.01 0.01 0.01 0.01 0.01
Cumulative Var 0.28 0.41 0.49 0.52 0.55 0.57 0.58 0.60 0.61 0.62 0.63 0.64
Proportion Explained 0.42 0.21 0.12 0.05 0.03 0.03 0.02 0.02 0.02 0.02 0.02 0.01
Cumulative Proportion 0.42 0.64 0.76 0.81 0.84 0.87 0.90 0.92 0.94 0.96 0.97 0.99
MR13
SS loadings 0.32
Proportion Var 0.01
Cumulative Var 0.65
Proportion Explained 0.01
Cumulative Proportion 1.00
Mean item complexity = 3.1
Test of the hypothesis that 13 factors are sufficient.
The degrees of freedom for the null model are 780 and the objective function was 27.45 with Chi Square of 5073.12
The degrees of freedom for the model are 338 and the objective function was 2.41
The root mean square of the residuals (RMSR) is 0.02
The df corrected root mean square of the residuals is 0.03
The harmonic number of observations is 196 with the empirical chi square 93.35 with prob < 1
The total number of observations was 200 with Likelihood Chi Square = 424.01 with prob < 0.001
Tucker Lewis Index of factoring reliability = 0.951
RMSEA index = 0.046 and the 90 % confidence intervals are 0.024 0.046
BIC = -1366.82
Fit based upon off diagonal values = 1
Measures of factor score adequacy
MR1 MR2 MR3 MR4 MR5 MR6 MR7
Correlation of scores with factors 0.99 0.98 0.95 0.88 0.89 0.84 0.83
Multiple R square of scores with factors 0.98 0.95 0.91 0.77 0.79 0.70 0.70
Minimum correlation of possible factor scores 0.95 0.90 0.82 0.55 0.59 0.40 0.39
MR8 MR9 MR10 MR11 MR12 MR13
Correlation of scores with factors 0.85 0.78 0.75 0.73 0.74 0.71
Multiple R square of scores with factors 0.73 0.61 0.57 0.54 0.54 0.50
Minimum correlation of possible factor scores 0.45 0.23 0.14 0.08 0.08 0.00
# count factors with eigenvalues > 1 and variance explained > 5%
efa_d1_all_unrotated_nfactors <- efa_d1_all_unrotated_eigenvalues %>%
filter(SS.loadings > 1, Proportion.Explained > 0.05) %>%
count() %>%
as.numeric()
efa_d1_all_unrotated_nfactors
[1] 3
efa_d1_rotatedN <- fa(d1_all, efa_d1_all_unrotated_nfactors, rotate = chosenRotType,
cor = chosenCorType, fm = "minres")
# check that each of these factors is the dominant factor for at least one mental capacity item
efa_d1_rotatedN_loadings <- fa.sort(loadings(efa_d1_rotatedN)[]) %>%
data.frame() %>%
rownames_to_column("capacity") %>%
gather(factor, loading, -capacity) %>%
mutate(loading_abs = abs(loading)) %>%
group_by(capacity) %>%
top_n(1, loading_abs) %>%
ungroup()
efa_d1_rotatedN_loadings
# drop any factors where n < 1
efa_d1_rotatedN_loadings %>%
count(factor) %>%
filter(n > 0)
# set number of factors to extract
nfactors_d1_all <- efa_d1_rotatedN_loadings %>%
count(factor) %>%
filter(n > 0) %>%
nrow()
nfactors_d1_all
[1] 3
# run EFA with rotation with N factors
efa_d1_rotatedN <- fa(d1_all, nfactors_d1_all,
rotate = chosenRotType, cor = chosenCorType, fm = "minres")
fa.sort(efa_d1_rotatedN)
Factor Analysis using method = minres
Call: fa(r = d1_all, nfactors = nfactors_d1_all, rotate = chosenRotType,
fm = "minres", cor = chosenCorType)
Standardized loadings (pattern matrix) based upon correlation matrix
MR1 MR2 MR3 h2 u2 com
pride 0.85 -0.06 -0.04 0.68 0.32 1.0
joy 0.84 0.12 -0.09 0.73 0.27 1.1
depressed 0.81 0.10 -0.06 0.68 0.32 1.0
happy 0.78 0.17 -0.04 0.68 0.32 1.1
love 0.76 0.10 0.06 0.66 0.34 1.0
guilt 0.75 -0.03 -0.05 0.53 0.47 1.0
disrespected 0.72 0.05 -0.06 0.53 0.47 1.0
embarrassed 0.66 -0.06 -0.07 0.39 0.61 1.0
beliefs 0.53 -0.17 0.24 0.40 0.60 1.6
personality 0.52 -0.18 0.19 0.34 0.66 1.5
angry 0.50 0.27 0.09 0.43 0.57 1.6
morality 0.43 -0.32 0.28 0.36 0.64 2.6
emo_recog 0.41 -0.18 0.26 0.29 0.71 2.1
self_restraint 0.40 -0.10 0.29 0.30 0.70 2.0
hungry 0.00 0.93 -0.06 0.87 0.13 1.0
computations 0.04 -0.83 0.33 0.80 0.20 1.3
pain 0.13 0.79 0.05 0.70 0.30 1.1
tired 0.17 0.70 0.20 0.64 0.36 1.3
fear 0.26 0.70 0.11 0.68 0.32 1.3
odors -0.08 0.69 0.25 0.50 0.50 1.3
safe 0.22 0.65 0.22 0.63 0.37 1.5
nauseated 0.23 0.62 0.13 0.54 0.46 1.4
desires 0.31 0.45 0.17 0.44 0.56 2.1
calm 0.39 0.45 0.12 0.50 0.50 2.1
pleasure 0.44 0.44 0.02 0.51 0.49 2.0
depth -0.17 0.15 0.62 0.35 0.65 1.3
reasoning 0.07 -0.02 0.60 0.39 0.61 1.0
seeing -0.17 0.29 0.59 0.37 0.63 1.6
choices 0.03 0.08 0.59 0.36 0.64 1.0
recognizing 0.14 -0.48 0.57 0.59 0.41 2.1
remembering 0.11 -0.37 0.56 0.48 0.52 1.8
temperature -0.13 0.21 0.55 0.30 0.70 1.4
conscious -0.05 0.34 0.54 0.37 0.63 1.7
communicating 0.07 -0.33 0.53 0.41 0.59 1.7
sounds -0.13 0.18 0.53 0.27 0.73 1.3
free_will 0.02 0.06 0.49 0.25 0.75 1.0
intentions 0.25 -0.38 0.43 0.41 0.59 2.6
goal 0.15 0.12 0.41 0.25 0.75 1.4
self_aware 0.15 0.18 0.40 0.27 0.73 1.7
thoughts 0.31 0.16 0.34 0.33 0.67 2.4
MR1 MR2 MR3
SS loadings 7.56 6.63 5.04
Proportion Var 0.19 0.17 0.13
Cumulative Var 0.19 0.35 0.48
Proportion Explained 0.39 0.34 0.26
Cumulative Proportion 0.39 0.74 1.00
With factor correlations of
MR1 MR2 MR3
MR1 1.00 0.28 0.30
MR2 0.28 1.00 -0.01
MR3 0.30 -0.01 1.00
Mean item complexity = 1.5
Test of the hypothesis that 3 factors are sufficient.
The degrees of freedom for the null model are 780 and the objective function was 27.45 with Chi Square of 5073.12
The degrees of freedom for the model are 663 and the objective function was 6.67
The root mean square of the residuals (RMSR) is 0.05
The df corrected root mean square of the residuals is 0.05
The harmonic number of observations is 196 with the empirical chi square 729.79 with prob < 0.036
The total number of observations was 200 with Likelihood Chi Square = 1219.34 with prob < 2.1e-35
Tucker Lewis Index of factoring reliability = 0.846
RMSEA index = 0.071 and the 90 % confidence intervals are 0.059 NA
BIC = -2293.44
Fit based upon off diagonal values = 0.98
Measures of factor score adequacy
MR1 MR2 MR3
Correlation of scores with factors 0.97 0.98 0.95
Multiple R square of scores with factors 0.95 0.96 0.89
Minimum correlation of possible factor scores 0.89 0.92 0.79
# get loadings for each factor
efa_d1_rotatedN_loadings <- loadings(efa_d1_rotatedN)[] %>%
data.frame() %>%
rownames_to_column(var = "capacity")
data.frame(loadings(fa.sort(efa_d1_rotatedN))[]) %>%
rownames_to_column("capacity") %>%
mutate_at(vars(starts_with("M")), funs(round2))
Study information:
# make demographics tables
demoSampleSize("study 2")
demoDuration("study 2")
demoAge("study 2")
Joining, by = c("character", "min_age", "max_age", "median_age", "mean_age", "sd_age")
Column `character` joining factor and character vector, coercing into character vector
demoGender("study 2")
demoRace("study 2") %>%
filter(character == "all") %>%
gather(ethnicity, count, -character) %>%
mutate(prop = round(count/sum(count), 2)) %>%
arrange(desc(prop))
# gender
with(d2 %>%
distinct(subid, character, gender) %>%
filter(!is.na(gender)) %>%
mutate(character = factor(character)),
table(character, gender)) %>%
summary()
Number of cases in table: 197
Number of factors: 2
Test for independence of all factors:
Chisq = 0.016177, df = 1, p-value = 0.8988
# age
with(d2 %>%
distinct(subid, character, age),
t.test(age ~ character))
Welch Two Sample t-test
data: age by character
t = 0.48254, df = 197.76, p-value = 0.63
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-0.1769700 0.2916335
sample estimates:
mean in group beetle mean in group robot
8.381464 8.324132
# examine scree plot
# fa.parallel(d2_all)
# run EFA without rotation with N factors
efa_d2_all_unrotated <- fa(d2_all, 13, rotate = "none",
cor = chosenCorType, fm = "minres")
fa.sort(efa_d2_all_unrotated)
Factor Analysis using method = minres
Call: fa(r = d2_all, nfactors = 13, rotate = "none", fm = "minres",
cor = chosenCorType)
Standardized loadings (pattern matrix) based upon correlation matrix
MR1 MR2 MR3 MR4 MR5 MR6 MR7 MR8 MR9 MR10 MR11
happy 0.72 0.09 -0.20 0.18 -0.20 -0.06 -0.01 -0.13 0.09 -0.03 0.12
joy 0.70 0.00 -0.26 0.10 -0.18 -0.10 -0.03 -0.05 0.06 -0.01 0.05
depressed 0.69 -0.11 -0.13 0.05 -0.06 -0.13 0.12 0.06 0.21 -0.09 -0.15
disrespected 0.69 0.00 -0.11 -0.07 -0.03 -0.02 0.16 -0.16 -0.09 0.03 0.10
pride 0.68 0.14 -0.31 0.06 -0.05 -0.02 0.01 -0.04 -0.01 -0.04 -0.03
love 0.60 0.03 -0.20 0.11 0.03 -0.11 0.09 -0.19 -0.05 0.05 0.07
pleasure 0.60 0.04 -0.18 0.10 -0.23 -0.17 -0.31 0.04 0.01 -0.04 0.08
safe 0.58 -0.08 0.23 0.06 -0.06 -0.06 -0.19 0.16 -0.24 0.00 -0.10
guilt 0.58 0.07 -0.18 0.13 0.25 0.19 -0.07 -0.10 -0.07 -0.08 0.14
thoughts 0.57 0.01 0.18 -0.05 -0.11 -0.16 -0.03 0.03 -0.09 0.03 -0.03
angry 0.57 -0.17 -0.04 -0.07 0.05 0.16 0.11 0.06 -0.04 0.01 -0.04
fear 0.55 -0.37 0.11 0.11 -0.08 0.03 0.01 0.17 0.09 0.03 -0.06
calm 0.55 -0.07 0.01 0.02 -0.10 -0.09 -0.14 0.05 -0.01 -0.18 0.10
personality 0.54 0.29 -0.01 -0.12 0.20 0.12 0.01 0.18 -0.01 0.06 -0.09
beliefs 0.54 0.23 -0.07 -0.09 0.05 0.24 -0.12 -0.05 -0.01 0.13 0.06
embarrassed 0.53 0.03 -0.30 0.10 0.22 0.34 0.23 -0.06 -0.01 -0.03 0.02
desires 0.52 -0.19 -0.03 -0.01 -0.08 0.21 -0.07 -0.04 0.02 0.10 0.03
free_will 0.49 0.00 0.32 -0.22 -0.09 0.11 -0.21 -0.18 0.08 0.03 -0.02
morality 0.44 0.37 -0.13 -0.08 0.21 -0.26 -0.05 0.26 -0.18 0.33 0.07
choices 0.43 0.05 0.26 -0.38 -0.07 0.20 -0.05 -0.13 0.11 0.05 -0.15
tired 0.39 -0.35 0.07 0.11 0.16 -0.16 0.24 0.06 0.12 -0.07 -0.11
goal 0.35 0.31 -0.03 -0.09 -0.06 0.10 -0.11 0.22 0.02 -0.18 -0.25
computations -0.01 0.80 -0.04 0.02 0.08 0.10 -0.09 0.03 0.10 -0.11 0.08
hungry 0.38 -0.77 0.22 0.02 0.14 -0.04 0.01 -0.08 0.06 0.05 0.00
pain 0.45 -0.64 0.21 -0.06 0.07 -0.09 -0.05 -0.04 -0.01 0.10 -0.01
remembering 0.06 0.58 0.17 0.13 -0.05 0.00 0.01 -0.17 0.19 0.28 -0.24
odors 0.15 -0.53 0.38 0.06 0.09 0.12 -0.06 0.16 0.06 0.03 0.19
emo_recog 0.34 0.50 -0.05 0.09 0.11 -0.07 0.22 0.14 -0.09 0.07 -0.02
nauseated 0.30 -0.43 0.06 0.11 0.27 0.05 0.00 0.07 0.16 0.06 -0.10
intentions 0.35 0.40 0.15 -0.24 0.11 0.03 0.02 -0.12 0.09 -0.14 0.08
recognizing 0.20 0.32 0.13 0.10 0.11 0.00 -0.11 0.21 0.06 0.00 0.04
communicating 0.09 0.30 0.16 0.30 -0.17 0.03 0.19 0.11 0.29 0.04 0.04
conscious 0.36 0.08 0.48 0.21 -0.14 0.01 0.16 -0.12 -0.27 -0.14 -0.04
self_aware 0.27 0.21 0.46 -0.03 -0.10 0.09 0.23 0.00 -0.31 -0.10 -0.09
temperature -0.05 0.34 0.41 0.31 0.34 -0.26 -0.15 -0.23 0.08 -0.15 -0.09
sounds -0.07 0.10 0.40 0.10 -0.18 0.23 0.08 0.23 0.14 -0.08 0.25
depth 0.14 0.22 0.36 0.12 0.34 0.02 -0.14 -0.07 -0.02 0.02 0.17
reasoning 0.23 0.28 0.36 0.05 -0.16 -0.01 -0.05 0.01 0.05 0.12 -0.06
self_restraint 0.34 0.19 0.15 -0.56 0.11 -0.33 0.20 0.05 0.17 -0.16 0.15
seeing -0.07 0.13 0.26 -0.03 -0.19 -0.08 0.19 -0.08 0.02 0.27 0.20
MR12 MR13 h2 u2 com
happy 0.05 -0.09 0.69 0.31 1.8
joy 0.11 0.02 0.62 0.38 1.6
depressed -0.01 0.00 0.63 0.37 1.7
disrespected 0.09 0.08 0.58 0.42 1.4
pride 0.01 -0.15 0.61 0.39 1.7
love 0.03 0.14 0.50 0.50 1.9
pleasure -0.09 -0.02 0.60 0.40 2.5
safe 0.09 -0.22 0.59 0.41 2.8
guilt -0.05 -0.03 0.53 0.47 2.4
thoughts 0.15 0.08 0.44 0.56 1.8
angry 0.02 -0.16 0.44 0.56 1.7
fear -0.25 -0.09 0.58 0.42 2.9
calm -0.02 0.10 0.40 0.60 1.8
personality 0.00 0.09 0.50 0.50 2.6
beliefs 0.08 0.21 0.50 0.50 2.8
embarrassed -0.13 -0.05 0.63 0.37 3.7
desires -0.17 0.15 0.42 0.58 2.3
free_will -0.03 -0.10 0.51 0.49 3.5
morality 0.03 -0.09 0.69 0.31 5.9
choices 0.16 -0.09 0.54 0.46 4.6
tired 0.12 0.05 0.45 0.55 4.7
goal -0.04 0.14 0.42 0.58 5.4
computations 0.10 0.00 0.71 0.29 1.2
hungry 0.03 0.01 0.81 0.19 1.8
pain 0.06 0.08 0.70 0.30 2.3
remembering -0.03 -0.08 0.59 0.41 2.8
odors 0.09 -0.06 0.55 0.45 3.0
emo_recog 0.13 0.05 0.49 0.51 3.1
nauseated 0.02 0.08 0.42 0.58 3.6
intentions 0.05 -0.12 0.43 0.57 4.2
recognizing -0.14 0.05 0.26 0.74 4.7
communicating -0.02 -0.12 0.40 0.60 5.9
conscious -0.12 0.01 0.58 0.42 4.2
self_aware -0.03 0.01 0.52 0.48 4.2
temperature 0.08 0.04 0.68 0.32 6.4
sounds 0.23 0.06 0.47 0.53 5.6
depth -0.15 -0.07 0.41 0.59 4.8
reasoning -0.11 0.22 0.37 0.63 4.8
self_restraint -0.19 0.01 0.76 0.24 4.5
seeing -0.08 0.02 0.30 0.70 6.0
MR1 MR2 MR3 MR4 MR5 MR6 MR7 MR8 MR9 MR10 MR11 MR12
SS loadings 8.25 4.13 2.18 1.05 0.96 0.84 0.73 0.66 0.62 0.54 0.52 0.44
Proportion Var 0.21 0.10 0.05 0.03 0.02 0.02 0.02 0.02 0.02 0.01 0.01 0.01
Cumulative Var 0.21 0.31 0.36 0.39 0.41 0.44 0.45 0.47 0.49 0.50 0.51 0.52
Proportion Explained 0.39 0.19 0.10 0.05 0.05 0.04 0.03 0.03 0.03 0.03 0.02 0.02
Cumulative Proportion 0.39 0.58 0.68 0.73 0.78 0.82 0.85 0.88 0.91 0.94 0.96 0.98
MR13
SS loadings 0.40
Proportion Var 0.01
Cumulative Var 0.53
Proportion Explained 0.02
Cumulative Proportion 1.00
Mean item complexity = 3.4
Test of the hypothesis that 13 factors are sufficient.
The degrees of freedom for the null model are 780 and the objective function was 17.08 with Chi Square of 3157.08
The degrees of freedom for the model are 338 and the objective function was 1.79
The root mean square of the residuals (RMSR) is 0.02
The df corrected root mean square of the residuals is 0.03
The harmonic number of observations is 198 with the empirical chi square 142.71 with prob < 1
The total number of observations was 200 with Likelihood Chi Square = 314.7 with prob < 0.81
Tucker Lewis Index of factoring reliability = 1.024
RMSEA index = 0.016 and the 90 % confidence intervals are 0 0.018
BIC = -1476.13
Fit based upon off diagonal values = 0.99
Measures of factor score adequacy
MR1 MR2 MR3 MR4 MR5 MR6 MR7
Correlation of scores with factors 0.98 0.96 0.91 0.86 0.83 0.83 0.79
Multiple R square of scores with factors 0.95 0.93 0.83 0.74 0.69 0.68 0.62
Minimum correlation of possible factor scores 0.90 0.86 0.67 0.48 0.38 0.36 0.23
MR8 MR9 MR10 MR11 MR12 MR13
Correlation of scores with factors 0.77 0.77 0.75 0.73 0.71 0.67
Multiple R square of scores with factors 0.60 0.59 0.56 0.53 0.50 0.45
Minimum correlation of possible factor scores 0.19 0.18 0.12 0.06 0.00 -0.10
# count factors with eigenvalues > 1 and variance explained > 5%
efa_d2_all_unrotated_nfactors <- efa_d2_all_unrotated_eigenvalues %>%
filter(SS.loadings > 1, Proportion.Explained > 0.05) %>%
count() %>%
as.numeric()
efa_d2_all_unrotated_nfactors
[1] 3
efa_d2_rotatedN <- fa(d2_all, efa_d2_all_unrotated_nfactors, rotate = chosenRotType,
cor = chosenCorType, fm = "minres")
# check that each of these factors is the dominant factor for at least one mental capacity item
efa_d2_rotatedN_loadings <- fa.sort(loadings(efa_d2_rotatedN)[]) %>%
data.frame() %>%
rownames_to_column("capacity") %>%
gather(factor, loading, -capacity) %>%
mutate(loading_abs = abs(loading)) %>%
group_by(capacity) %>%
top_n(1, loading_abs) %>%
ungroup()
efa_d2_rotatedN_loadings
# drop any factors where n < 1
efa_d2_rotatedN_loadings %>%
count(factor) %>%
filter(n > 0)
# set number of factors to extract
nfactors_d2_all <- efa_d2_rotatedN_loadings %>%
count(factor) %>%
filter(n > 0) %>%
nrow()
nfactors_d2_all
[1] 3
# run EFA with rotation with N factors
efa_d2_rotatedN <- fa(d2_all, nfactors_d2_all,
rotate = chosenRotType, cor = chosenCorType, fm = "minres")
fa.sort(efa_d2_rotatedN)
Factor Analysis using method = minres
Call: fa(r = d2_all, nfactors = nfactors_d2_all, rotate = chosenRotType,
fm = "minres", cor = chosenCorType)
Standardized loadings (pattern matrix) based upon correlation matrix
MR1 MR2 MR3 h2 u2 com
pride 0.81 -0.11 -0.08 0.594 0.41 1.1
joy 0.76 0.05 -0.07 0.555 0.44 1.0
happy 0.74 0.00 0.02 0.555 0.44 1.0
disrespected 0.65 0.11 0.08 0.488 0.51 1.1
depressed 0.65 0.20 0.02 0.504 0.50 1.2
love 0.64 0.03 -0.03 0.403 0.60 1.0
embarrassed 0.62 -0.02 -0.12 0.346 0.65 1.1
pleasure 0.61 0.03 0.00 0.375 0.63 1.0
guilt 0.60 -0.01 0.01 0.358 0.64 1.0
beliefs 0.51 -0.12 0.16 0.336 0.66 1.3
personality 0.49 -0.15 0.24 0.372 0.63 1.7
angry 0.48 0.26 0.05 0.356 0.64 1.6
morality 0.47 -0.26 0.13 0.306 0.69 1.8
calm 0.44 0.17 0.12 0.299 0.70 1.5
desires 0.42 0.27 0.04 0.301 0.70 1.7
thoughts 0.35 0.17 0.34 0.364 0.64 2.5
goal 0.35 -0.21 0.17 0.212 0.79 2.2
hungry 0.06 0.87 0.04 0.786 0.21 1.0
pain 0.15 0.77 0.10 0.662 0.34 1.1
computations 0.12 -0.76 0.25 0.650 0.35 1.3
odors -0.20 0.64 0.22 0.425 0.57 1.4
fear 0.33 0.49 0.12 0.441 0.56 1.9
nauseated 0.15 0.47 -0.02 0.270 0.73 1.2
remembering 0.02 -0.44 0.36 0.333 0.67 1.9
tired 0.23 0.42 0.04 0.266 0.73 1.6
emo_recog 0.38 -0.40 0.22 0.360 0.64 2.6
self_aware -0.05 0.02 0.57 0.303 0.70 1.0
conscious 0.00 0.16 0.56 0.331 0.67 1.2
reasoning -0.02 -0.08 0.51 0.265 0.73 1.0
depth -0.09 -0.04 0.45 0.182 0.82 1.1
free_will 0.19 0.21 0.43 0.331 0.67 1.9
temperature -0.23 -0.17 0.43 0.212 0.79 1.9
sounds -0.28 0.03 0.37 0.144 0.86 1.9
intentions 0.24 -0.24 0.37 0.293 0.71 2.5
choices 0.19 0.13 0.36 0.231 0.77 1.8
safe 0.31 0.27 0.34 0.380 0.62 2.9
recognizing 0.12 -0.20 0.29 0.154 0.85 2.1
seeing -0.21 -0.04 0.27 0.083 0.92 1.9
communicating 0.01 -0.20 0.26 0.111 0.89 1.9
self_restraint 0.21 -0.04 0.26 0.144 0.86 2.0
MR1 MR2 MR3
SS loadings 6.84 4.15 3.09
Proportion Var 0.17 0.10 0.08
Cumulative Var 0.17 0.27 0.35
Proportion Explained 0.49 0.29 0.22
Cumulative Proportion 0.49 0.78 1.00
With factor correlations of
MR1 MR2 MR3
MR1 1.00 0.16 0.34
MR2 0.16 1.00 -0.02
MR3 0.34 -0.02 1.00
Mean item complexity = 1.6
Test of the hypothesis that 3 factors are sufficient.
The degrees of freedom for the null model are 780 and the objective function was 17.08 with Chi Square of 3157.08
The degrees of freedom for the model are 663 and the objective function was 4.79
The root mean square of the residuals (RMSR) is 0.05
The df corrected root mean square of the residuals is 0.05
The harmonic number of observations is 198 with the empirical chi square 788.89 with prob < 0.00052
The total number of observations was 200 with Likelihood Chi Square = 874.92 with prob < 5.4e-08
Tucker Lewis Index of factoring reliability = 0.894
RMSEA index = 0.047 and the 90 % confidence intervals are 0.032 0.047
BIC = -2637.86
Fit based upon off diagonal values = 0.95
Measures of factor score adequacy
MR1 MR2 MR3
Correlation of scores with factors 0.96 0.96 0.90
Multiple R square of scores with factors 0.92 0.92 0.81
Minimum correlation of possible factor scores 0.85 0.84 0.62
# get loadings for each factor
efa_d2_rotatedN_loadings <- loadings(efa_d2_rotatedN)[] %>%
data.frame() %>%
rownames_to_column(var = "capacity")
data.frame(loadings(fa.sort(efa_d2_rotatedN))[]) %>%
rownames_to_column("capacity") %>%
mutate_at(vars(starts_with("M")), funs(round2))
# manually set 3 factors
order_s1_manual <- loadings(fa.sort(fa(d1_all, nfactors = 3,
rotate = chosenRotType, cor = chosenCorType)))[] %>%
data.frame() %>%
rownames_to_column(var = "capacity") %>%
rownames_to_column(var = "order1_manual") %>%
rename(s1_heart = MR1, s1_body = MR2, s1_mind = MR3)
order_s1 <- loadings(fa.sort(efa_d1_rotatedN))[] %>%
data.frame() %>%
rownames_to_column(var = "capacity") %>%
rownames_to_column(var = "order1") %>%
rename(s1_heart = MR1, s1_body = MR2, s1_mind = MR3)
order_s2 <- loadings(fa.sort(efa_d2_rotatedN))[] %>%
data.frame() %>%
rownames_to_column(var = "capacity") %>%
rename(s2_body = MR2, s2_heart = MR1, s2_mind = MR3)
bigTable <- order_s1 %>%
# order_s1_manual %>% # could substitute
full_join(order_s2) %>%
mutate_at(vars(starts_with("s")), funs(round2)) %>%
select(order1, # could subistitute order_s1
# select(order1_manual, # could subistitute order_s1
capacity, ends_with("heart"), ends_with("body"), ends_with("mind"))
Joining, by = "capacity"
bigTable
Factor loadings for the 40 mental capacities on the three rotated factors in Study 1. Items are colored by their dominant factor loading: Items that loaded most strongly on the body factor (bodily states and will) are in red; items that loaded most strongly on the heart factor (social-emotional experiences and morality) are in blue; and items that loaded most strongly on the mind factor (perceptual-cognitive abilities and goal pursuit) are in green.
# set up labels for plot (shortened version of mental capacity items)
wording_s1 <- loadings(efa_d1_rotatedN)[] %>%
data.frame() %>%
rownames_to_column(var = "item") %>%
select(item) %>%
mutate(wording = factor(
recode(item,
happy = "feel happy",
depressed = "feel sad",
fear = "feel scared",
angry = "get angry",
calm = "feel calm",
sounds = "hear sounds",
seeing = "see things",
temperature = "sense temperatures",
odors = "smell things",
depth = "sense whether something is close by or far away",
computations = "do math",
thoughts = "have thoughts",
reasoning = "figure out how to do things",
remembering = "remember things",
beliefs = "have beliefs, like when you think something is true",
hungry = "get hungry",
tired = "feel tired",
pain = "feel pain",
nauseated = "feel sick, like when you feel like you might throw up",
safe = "feel safe",
love = "feel love",
recognizing = "recognize somebody else",
communicating = "communicate with somebody else",
guilt = "feel guilty",
disrespected = "get hurt feelings",
free_will = "decide what to do",
choices = "make choices",
self_restraint = "have self-control, like when you stop yourself from doing something you shouldn't do",
intentions = "make plans",
goal = "have goals, like when you're working hard to do something or make something happen",
conscious = "be aware of things",
self_aware = "be aware of itself",
desires = "have desires, like when you really want something",
embarrassed = "feel embarrassed",
emo_recog = "understand how somebody else is feeling",
joy = "feel joy",
morality = "know what's nice and what's mean",
personality = "have a personality, like when someone is shy and somebody else is silly",
pleasure = "feel pleasure, like when something feels really good",
pride = "feel proud"))) %>%
mutate(short = factor(
recode(item,
happy = "happy",
depressed = "sad",
fear = "scared",
angry = "angry",
calm = "calm",
sounds = "hear",
seeing = "see",
temperature = "temperatures",
odors = "smell",
depth = "depth",
computations = "math",
thoughts = "thoughts",
reasoning = "figure out",
remembering = "remember",
beliefs = "beliefs",
hungry = "hungry",
tired = "tired",
pain = "pain",
nauseated = "sick",
safe = "safe",
love = "love",
recognizing = "recognize",
communicating = "communicate",
guilt = "guilty",
disrespected = "hurt feelings",
free_will = "decide",
choices = "choices",
self_restraint = "self-control",
intentions = "plans",
goal = "goals",
conscious = "aware",
self_aware = "self-aware",
desires = "desires",
embarrassed = "embarrassed",
emo_recog = "empathy",
joy = "joy",
morality = "morality",
personality = "personality",
pleasure = "pleasure",
pride = "pride")))
# make dataframe for plotting
scatter_plotting <- loadings(efa_d1_rotatedN)[] %>%
data.frame() %>%
rownames_to_column(var = "item") %>%
rename(HEART = MR1,
BODY = MR2,
MIND = MR3) %>%
full_join(wording_s1) %>%
mutate(dominant = factor(
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(BODY), "BODY",
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(HEART), "HEART",
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(MIND), "MIND",
NA)))),
size = ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(BODY), abs(BODY),
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(HEART), abs(HEART),
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(MIND), abs(MIND),
NA))),
color = ifelse(dominant == "BODY", "#E41A1C",
ifelse(dominant == "HEART", "#377EB8",
ifelse(dominant == "MIND", "#4DAF4A",
NA))))
Joining, by = "item"
# plot!
figS1 <- plot_ly(scatter_plotting, x = ~HEART, y = ~BODY, z = ~MIND,
type = "scatter3d",
color = ~dominant, colors = c("#E41A1C", "#377EB8", "#4DAF4A"),
marker = list(size = 4),
text = ~short,
textfont = list(size = 15),
mode = "text+markers",
showlegend = TRUE)
figS1
# set up labels for plot (shortened version of mental capacity items)
wording_s2 <- loadings(efa_d2_rotatedN)[] %>%
data.frame() %>%
rownames_to_column(var = "item") %>%
select(item) %>%
mutate(wording = factor(
recode(item,
happy = "feel happy",
depressed = "feel sad",
fear = "feel scared",
angry = "get angry",
calm = "feel calm",
sounds = "hear sounds",
seeing = "see things",
temperature = "sense temperatures",
odors = "smell things",
depth = "sense whether something is close by or far away",
computations = "do math",
thoughts = "have thoughts",
reasoning = "figure out how to do things",
remembering = "remember things",
beliefs = "have beliefs, like when you think something is true",
hungry = "get hungry",
tired = "feel tired",
pain = "feel pain",
nauseated = "feel sick, like when you feel like you might throw up",
safe = "feel safe",
love = "feel love",
recognizing = "recognize somebody else",
communicating = "communicate with somebody else",
guilt = "feel guilty",
disrespected = "get hurt feelings",
free_will = "decide what to do",
choices = "make choices",
self_restraint = "have self-control, like when you stop yourself from doing something you shouldn't do",
intentions = "make plans",
goal = "have goals, like when you're working hard to do something or make something happen",
conscious = "be aware of things",
self_aware = "be aware of itself",
desires = "have desires, like when you really want something",
embarrassed = "feel embarrassed",
emo_recog = "understand how somebody else is feeling",
joy = "feel joy",
morality = "know what's nice and what's mean",
personality = "have a personality, like when someone is shy and somebody else is silly",
pleasure = "feel pleasure, like when something feels really good",
pride = "feel proud"))) %>%
mutate(short = factor(
recode(item,
happy = "happy",
depressed = "sad",
fear = "scared",
angry = "angry",
calm = "calm",
sounds = "hear",
seeing = "see",
temperature = "temperatures",
odors = "smell",
depth = "depth",
computations = "math",
thoughts = "thoughts",
reasoning = "figure out",
remembering = "remember",
beliefs = "beliefs",
hungry = "hungry",
tired = "tired",
pain = "pain",
nauseated = "sick",
safe = "safe",
love = "love",
recognizing = "recognize",
communicating = "communicate",
guilt = "guilty",
disrespected = "hurt feelings",
free_will = "decide",
choices = "choices",
self_restraint = "self-control",
intentions = "plans",
goal = "goals",
conscious = "aware",
self_aware = "self-aware",
desires = "desires",
embarrassed = "embarrassed",
emo_recog = "empathy",
joy = "joy",
morality = "morality",
personality = "personality",
pleasure = "pleasure",
pride = "pride")))
# make dataframe for plotting
scatter_plotting <- loadings(efa_d2_rotatedN)[] %>%
data.frame() %>%
rownames_to_column(var = "item") %>%
rename(HEART = MR1,
BODY = MR2,
MIND = MR3) %>%
full_join(wording_s2) %>%
mutate(dominant = factor(
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(BODY), "BODY",
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(HEART), "HEART",
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(MIND), "MIND",
NA)))),
size = ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(BODY), abs(BODY),
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(HEART), abs(HEART),
ifelse(pmax(abs(BODY), abs(HEART), abs(MIND)) == abs(MIND), abs(MIND),
NA))),
color = ifelse(dominant == "BODY", "#E41A1C",
ifelse(dominant == "HEART", "#4DAF4A",
ifelse(dominant == "MIND", "#E41A1C",
NA))))
Joining, by = "item"
# plot!
figS2 <- plot_ly(scatter_plotting, x = ~HEART, y = ~BODY, z = ~MIND,
type = "scatter3d",
color = ~dominant, colors = c("#E41A1C", "#377EB8", "#4DAF4A"),
marker = list(size = 4),
text = ~short,
textfont = list(size = 15),
mode = "text+markers",
showlegend = TRUE)
figS2
NOTE: set to 3 factors manually, for now.
factors_s1 <- fa.sort(fa(d1_all, nfactors = 3, cor = chosenCorType, rotate = chosenRotType)$loadings[]) %>%
data.frame() %>%
rownames_to_column(var = "item") %>%
full_join(wording_s1) %>%
select(wording, MR1, MR2, MR3) %>%
rename(capacity = wording, Factor1 = MR1, Factor2 = MR2, Factor3 = MR3) %>%
rownames_to_column(var = "order") %>%
mutate(order = as.numeric(order))
Joining, by = "item"
factors_s1_long <- factors_s1 %>%
gather(factor, loading, -capacity, -order) %>%
mutate(factor = factor(gsub("Factor", "F", factor))) %>%
# mutate(factor = factor(gsub("Factor", "F", factor),
# levels = c("F1", "F3", "F2"))) %>%
# mutate(factor = factor(gsub("Factor", "F", factor),
# levels = c("F2", "F1", "F3"))) %>%
arrange(order, factor)
factors_s1_blank1 <- factors_s1_long %>%
mutate(loading = rep(100, length(factors_s1_long$loading)))
# factors_s1_blank2 <- factors_s1_long %>%
# mutate(loading = ifelse(factor == "F1", loading, rep(100, length(factors_s1_long$loading)*2/3)))
factors_s1_blank2 <- factors_s1_long %>%
mutate(loading = ifelse(factor == "F2", loading, rep(100, length(factors_s1_long$loading)*2/3)))
factors_s1_blank3 <- factors_s1_long %>%
mutate(loading = ifelse(factor != "F3", loading, rep(100, length(factors_s1_long$loading)*1/3)))
# ggplot(factors_s1_blank1, aes(x = factor,
# y = reorder(capacity, desc(order)), fill = loading)) +
# geom_tile(color = "black") +
# # geom_text(aes(label = format(round(loading, 2), nsmall = 2))) +
# scale_fill_distiller(palette = "RdYlBu", limits = c(-1, 1), breaks = c(-1, 0, 1),
# guide = guide_colorbar(title = element_blank(),
# barheight = 20)) +
# scale_x_discrete(position = "top") +
# theme_minimal() +
# theme(text = element_text(size = 24),
# axis.text.x = element_text(size = 28),
# axis.title = element_blank(),
# panel.grid = element_blank()) # 1000 by 1000
#
# ggplot(factors_s1_blank2, aes(x = factor,
# y = reorder(capacity, desc(order)), fill = loading)) +
# geom_tile(color = "black") +
# # geom_text(aes(label = format(round(loading, 2), nsmall = 2))) +
# scale_fill_distiller(palette = "RdYlBu", limits = c(-1, 1), breaks = c(-1, 0, 1),
# guide = guide_colorbar(title = element_blank(),
# barheight = 20)) +
# scale_x_discrete(position = "top") +
# theme_minimal() +
# theme(text = element_text(size = 24),
# axis.text.x = element_text(size = 28),
# axis.title = element_blank(),
# panel.grid = element_blank()) # 1000 by 1000
#
# ggplot(factors_s1_blank3, aes(x = factor,
# y = reorder(capacity, desc(order)), fill = loading)) +
# geom_tile(color = "black") +
# # geom_text(aes(label = format(round(loading, 2), nsmall = 2))) +
# scale_fill_distiller(palette = "RdYlBu", limits = c(-1, 1), breaks = c(-1, 0, 1),
# guide = guide_colorbar(title = element_blank(),
# barheight = 20)) +
# scale_x_discrete(position = "top") +
# theme_minimal() +
# theme(text = element_text(size = 24),
# axis.text.x = element_text(size = 28),
# axis.title = element_blank(),
# panel.grid = element_blank()) # 1000 by 1000
ggplot(factors_s1_long, aes(x = factor,
y = reorder(capacity, desc(order)), fill = loading)) +
geom_tile(color = "black") +
geom_text(aes(label = format(round(loading, 2), nsmall = 2)), size = 6) +
scale_fill_distiller(palette = "RdYlBu", limits = c(-1, 1), breaks = c(-1, 0, 1),
guide = guide_colorbar(title = element_blank(),
barheight = 20)) +
scale_x_discrete(position = "top") +
# geom_rect(aes(xmin = 0.51, xmax = 1.49, ymin = 14.55, ymax = 20.45),
# alpha = 0, color = "black", size = .5) +
# geom_rect(aes(xmin = 1.51, xmax = 2.49, ymin = 6.55, ymax = 14.45),
# alpha = 0, color = "black", size = .5) +
# geom_rect(aes(xmin = 2.51, xmax = 3.49, ymin = 0.55, ymax = 6.45),
# alpha = 0, color = "black", size = .5) +
# theme_bw() +
theme_minimal() +
theme(text = element_text(size = 24),
axis.text.x = element_text(size = 28),
axis.title = element_blank(),
panel.grid = element_blank()) # 1000 by 1000
factors_s2 <- fa.sort(fa(d2_all, nfactors = 3, cor = chosenCorType, rotate = chosenRotType)$loadings[]) %>%
data.frame() %>%
rownames_to_column(var = "item") %>%
full_join(wording_s2) %>%
select(wording, MR1, MR2, MR3) %>%
rename(capacity = wording, Factor1 = MR1, Factor2 = MR2, Factor3 = MR3) %>%
rownames_to_column(var = "order") %>%
mutate(order = as.numeric(order))
Joining, by = "item"
factors_s2_long <- factors_s2 %>%
gather(factor, loading, -capacity, -order) %>%
mutate(factor = factor(gsub("Factor", "F", factor))) %>%
# mutate(factor = factor(gsub("Factor", "F", factor),
# levels = c("F1", "F3", "F2"))) %>%
# mutate(factor = factor(gsub("Factor", "F", factor),
# levels = c("F2", "F1", "F3"))) %>%
arrange(order, factor)
factors_s2_blank1 <- factors_s2_long %>%
mutate(loading = rep(100, length(factors_s2_long$loading)))
factors_s2_blank2 <- factors_s2_long %>%
mutate(loading = ifelse(factor == "F1", loading, rep(100, length(factors_s2_long$loading)*2/3)))
# factors_s2_blank2 <- factors_s2_long %>%
# mutate(loading = ifelse(factor == "F2", loading, rep(100, length(factors_s2_long$loading)*2/3)))
factors_s2_blank3 <- factors_s2_long %>%
mutate(loading = ifelse(factor != "F3", loading, rep(100, length(factors_s2_long$loading)*1/3)))
# ggplot(factors_s2_blank1, aes(x = factor,
# y = reorder(capacity, desc(order)), fill = loading)) +
# geom_tile(color = "black") +
# # geom_text(aes(label = format(round(loading, 2), nsmall = 2))) +
# scale_fill_distiller(palette = "RdYlBu", limits = c(-1, 1), breaks = c(-1, 0, 1),
# guide = guide_colorbar(title = element_blank(),
# barheight = 20)) +
# scale_x_discrete(position = "top") +
# theme_minimal() +
# theme(text = element_text(size = 24),
# axis.text.x = element_text(size = 28),
# axis.title = element_blank(),
# panel.grid = element_blank()) # 1000 by 1000
#
# ggplot(factors_s2_blank2, aes(x = factor,
# y = reorder(capacity, desc(order)), fill = loading)) +
# geom_tile(color = "black") +
# # geom_text(aes(label = format(round(loading, 2), nsmall = 2))) +
# scale_fill_distiller(palette = "RdYlBu", limits = c(-1, 1), breaks = c(-1, 0, 1),
# guide = guide_colorbar(title = element_blank(),
# barheight = 20)) +
# scale_x_discrete(position = "top") +
# theme_minimal() +
# theme(text = element_text(size = 24),
# axis.text.x = element_text(size = 28),
# axis.title = element_blank(),
# panel.grid = element_blank()) # 1000 by 1000
#
# ggplot(factors_s2_blank3, aes(x = factor,
# y = reorder(capacity, desc(order)), fill = loading)) +
# geom_tile(color = "black") +
# # geom_text(aes(label = format(round(loading, 2), nsmall = 2))) +
# scale_fill_distiller(palette = "RdYlBu", limits = c(-1, 1), breaks = c(-1, 0, 1),
# guide = guide_colorbar(title = element_blank(),
# barheight = 20)) +
# scale_x_discrete(position = "top") +
# theme_minimal() +
# theme(text = element_text(size = 24),
# axis.text.x = element_text(size = 28),
# axis.title = element_blank(),
# panel.grid = element_blank()) # 1000 by 1000
ggplot(factors_s2_long, aes(x = factor,
y = reorder(capacity, desc(order)), fill = loading)) +
geom_tile(color = "black") +
geom_text(aes(label = format(round(loading, 2), nsmall = 2)), size = 6) +
scale_fill_distiller(palette = "RdYlBu", limits = c(-1, 1), breaks = c(-1, 0, 1),
guide = guide_colorbar(title = element_blank(),
barheight = 20)) +
scale_x_discrete(position = "top") +
# geom_rect(aes(xmin = 0.51, xmax = 1.49, ymin = 14.55, ymax = 20.45),
# alpha = 0, color = "black", size = .5) +
# geom_rect(aes(xmin = 1.51, xmax = 2.49, ymin = 6.55, ymax = 14.45),
# alpha = 0, color = "black", size = .5) +
# geom_rect(aes(xmin = 2.51, xmax = 3.49, ymin = 0.55, ymax = 6.45),
# alpha = 0, color = "black", size = .5) +
# theme_bw() +
theme_minimal() +
theme(text = element_text(size = 24),
axis.text.x = element_text(size = 28),
axis.title = element_blank(),
panel.grid = element_blank()) # 1000 by 1000
Mean ratings of 40 mental capacities for the 2 entities included in Studies 1-2. Participants responded on a 3-point scale (0 = “no”, 0.5 = “kinda”, 1 = “yes”). Error bars are nonparametric bootstrapped 95% confidence intervals. Mental capacities are grouped according to their dominant factor loading in Study 1 (adults).
# make dataframe
s12_plotting <- char_plotting %>%
filter(study %in% c("study 1", "study 2")) %>%
distinct()
# plot! (ordered by study 1 factor loadings)
s12 <- ggplot(s12_plotting,
aes(y = Mean, x = reorder(wording, desc(s1_order)),
colour = factor(s1_color), shape = study)) +
geom_point(stat = "identity", position = position_dodge(width = 0.6), size = 2) +
geom_errorbar(aes(ymin = Lower, ymax = Upper), width = 0.4,
position = position_dodge(width = 0.6)) +
facet_wrap(~ character) +
theme_bw() +
scale_y_continuous(name = "\nMean rating",
limits = c(0, 1),
breaks = c(0, 0.5, 1),
labels = c("0\n(no)", "0.5\n(kinda)", "1\n(yes)")) +
scale_shape_discrete(name = "Study:",
labels = c("Study 1: adults", "Study 2: 7-9y")) +
# scale_colour_brewer(name = "Factor:",
# type = "qual", palette = 6,
# guide = FALSE) +
scale_colour_manual(name = "Factor:",
values = c("#E41A1C", "#4DAF4A", "#377EB8"),
labels = c("BODY", "MIND", "HEART")) +
coord_flip() +
theme(text = element_text(size = 9),
axis.title.y = element_blank(),
axis.text.y = element_text(face = "italic",
colour = palette_s1),
panel.grid.minor = element_blank(),
legend.position = "right")
s12
r1b <- brm(score ~ character * factor * age_group + (1 | subid) , tempC,
family = "gaussian", cores = n_cores)
Compiling the C++ model
Start sampling
starting worker pid=1176 on localhost:11418 at 13:16:46.683
starting worker pid=1184 on localhost:11418 at 13:16:46.951
starting worker pid=1192 on localhost:11418 at 13:16:47.198
starting worker pid=1200 on localhost:11418 at 13:16:47.439
SAMPLING FOR MODEL 'gaussian(identity) brms-model' NOW (CHAIN 1).
Gradient evaluation took 0.012079 seconds
1000 transitions using 10 leapfrog steps per transition would take 120.79 seconds.
Adjust your expectations accordingly!
Iteration: 1 / 2000 [ 0%] (Warmup)
SAMPLING FOR MODEL 'gaussian(identity) brms-model' NOW (CHAIN 2).
Gradient evaluation took 0.015809 seconds
1000 transitions using 10 leapfrog steps per transition would take 158.09 seconds.
Adjust your expectations accordingly!
Iteration: 1 / 2000 [ 0%] (Warmup)
SAMPLING FOR MODEL 'gaussian(identity) brms-model' NOW (CHAIN 3).
Gradient evaluation took 0.022298 seconds
1000 transitions using 10 leapfrog steps per transition would take 222.98 seconds.
Adjust your expectations accordingly!
Iteration: 1 / 2000 [ 0%] (Warmup)
SAMPLING FOR MODEL 'gaussian(identity) brms-model' NOW (CHAIN 4).
Gradient evaluation took 0.025095 seconds
1000 transitions using 10 leapfrog steps per transition would take 250.95 seconds.
Adjust your expectations accordingly!
Iteration: 1 / 2000 [ 0%] (Warmup)
# r1c <- brm(score ~ character * factor * age_group + (factor | subid) , tempC,
# family = "gaussian", cores = n_cores)
# summary(r1c)
# plot
ggplot(scores_s12_plotting %>%
ungroup() %>%
mutate(factor = factor(factor,
# labels = c("Social-emotional",
# "Bodily",
# "Perceptual-cognitive")),
labels = c("Heart",
"Body",
"Mind")),
age_group = factor(age_group,
# levels = c("adults", "children_79"),
# labels = c("adults", "children"))),
levels = c("children_79", "adults"),
labels = c("children", "adults"))),
aes(x = age_group, y = Mean, color = character, shape = character)) +
facet_wrap("factor", ncol = 3) +
theme_bw() +
theme(text = element_text(size = 28),
legend.position = "bottom") +
geom_point(size = 5, position = position_dodge(width = 0.4)) +
geom_errorbar(aes(ymin = Lower, ymax = Upper),
width = 0, position = position_dodge(width = 0.4)) +
scale_shape_manual(values = c(19, 15)) +
labs(title = "Factor scores by age group",
# subtitle = "Adults (Study 1) vs. children (Study 2)\n",
x = "Age group",
y = "Mean factor score") # 1000 by 500
scores_s1_plotting <- d1 %>%
select(subid, age, character) %>%
distinct() %>%
mutate(subid = paste(character, subid, sep = "_")) %>%
full_join(efa_d1_rotatedN$scores %>%
data.frame() %>%
rownames_to_column("subid")) %>%
mutate(character = factor(character)) %>%
rename(score_F1 = MR1, score_F2 = MR2, score_F3 = MR3) %>%
filter(!is.na(score_F1), !is.na(score_F2), !is.na(score_F3), !is.na(age)) %>%
gather(factor, score, starts_with("score_")) %>%
mutate(factor = factor(factor))
ggplot(scores_s1_plotting %>%
ungroup() %>%
mutate(factor = factor(factor,
levels = c("score_F1", "score_F2", "score_F3"),
# labels = c("Social-emotional",
# "Bodily",
# "Perceptual-cognitive"))),
labels = c("Heart",
"Body",
"Mind"))),
aes(x = age, y = score, color = character, fill = character, shape = character)) +
facet_wrap("factor", ncol = 3) +
theme_bw() +
theme(text = element_text(size = 28),
legend.position = "bottom") +
# geom_smooth(method = "loess", alpha = 0.4) +
geom_smooth(method = "lm", alpha = 0.4) +
geom_point(size = 2) +
scale_shape_manual(values = c(19, 15)) +
labs(title = "Factor scores by adults' age",
# subtitle = "Adults (Study 1)\n",
x = "Age (years)",
y = "Factor score") # 1000 by 500
scores_s2_plotting <- d2 %>%
select(subid, age, character) %>%
distinct() %>%
mutate(subid = paste(character, subid, sep = "_")) %>%
# full_join(efa_d2_rotatedN$scores %>%
full_join(efa_d12_all_rotated$scores %>%
data.frame() %>%
rownames_to_column("subid")) %>%
mutate(character = factor(character)) %>%
rename(score_F1 = MR1, score_F2 = MR2, score_F3 = MR3) %>%
filter(!is.na(score_F1), !is.na(score_F2), !is.na(score_F3), !is.na(age)) %>%
gather(factor, score, starts_with("score_")) %>%
mutate(factor = factor(factor))
contrasts(scores_s2_plotting$factor) = cbind(factor1 = c(1, -1, 0),
factor3 = c(0, -1, 1))
contrasts(scores_s2_plotting$character) = cbind(robot = c(-1, 1))
r2 <- lmer(score ~ character * factor * scale(age, scale = F) + (1 | subid) , scores_s2_plotting)
summary(r2)
# r2b <- brm(score ~ character * factor * age + (1 | subid) , scores_s2_plotting %>% mutate(age = scale(age, scale = F)), family = "gaussian", cores = n_cores)
# summary(r2b)
# look at heart only
# r2c_heart <- brm(score ~ character * age,
# data = scores_s2_plotting %>%
# filter(factor == "score_F1") %>%
# mutate(age = scale(age, scale = F)),
# family = "gaussian", cores = n_cores)
# summary(r2c_heart)
# look at body only
# r2c_body <- brm(score ~ character * age,
# data = scores_s2_plotting %>%
# filter(factor == "score_F2") %>%
# mutate(age = scale(age, scale = F)),
# family = "gaussian", cores = n_cores)
# summary(r2c_body)
# look at mind only
# r2c_mind <- brm(score ~ character * age,
# data = scores_s2_plotting %>%
# filter(factor == "score_F3") %>%
# mutate(age = scale(age, scale = F)),
# family = "gaussian", cores = n_cores)
# summary(r2c_mind)
ggplot(scores_s2_plotting %>%
ungroup() %>%
mutate(factor = factor(factor,
labels = c("Heart",
"Body",
"Mind"))),
# labels = c("Social-emotional",
# "Bodily",
# "Perceptual-cognitive"))),
aes(x = age, y = score, color = character, fill = character, shape = character)) +
facet_wrap("factor", ncol = 3) +
theme_bw() +
theme(text = element_text(size = 28),
legend.position = "bottom") +
# geom_smooth(method = "loess", alpha = 0.4) +
geom_smooth(method = "lm", alpha = 0.4) +
geom_point(size = 2) +
scale_shape_manual(values = c(19, 15)) +
labs(title = "Factor scores by children's age",
# subtitle = "Children (Study 2)\n",
x = "Age (years)",
y = "Factor score") # 1000 by 500
tempA <- d2 %>%
select(subid, age, character) %>%
distinct() %>%
mutate(subid = paste(character, subid, sep = "_")) %>%
full_join(efa_d12_all_rotatedN$scores %>%
data.frame() %>%
rownames_to_column("subid")) %>%
mutate(character = factor(character)) %>%
rename(score_F1 = MR1, score_F2 = MR2, score_F3 = MR3) %>%
filter(!is.na(score_F1), !is.na(score_F2), !is.na(score_F3), !is.na(age)) %>%
gather(factor, score, starts_with("score_")) %>%
mutate(factor = factor(factor,
labels = c("Heart",
"Body",
"Mind")))
# labels = c("Social-emotional",
# "Bodily",
# "Perceptual-cognitive")))
tempB <- scores_s12_plotting %>%
filter(age_group == "adults") %>%
ungroup() %>%
mutate(factor = factor(factor,
labels = c("Heart",
"Body",
"Mind")),
# labels = c("Social-emotional",
# "Bodily",
# "Perceptual-cognitive")),
age = 11)
ggplot(tempA,
aes(x = age, y = score, color = character, fill = character, shape = character)) +
facet_wrap("factor", ncol = 3) +
theme_bw() +
theme(text = element_text(size = 28),
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1),
legend.position = "bottom") +
geom_hline(data = tempB, aes(yintercept = Mean, color = character), lty = 2) +
# geom_smooth(method = "loess", alpha = 0.4) +
geom_smooth(method = "lm", alpha = 0.4) +
geom_point(size = 2) +
geom_point(data = tempB, aes(y = Mean),
size = 4, position = position_dodge(width = 0.6)) +
geom_errorbar(data = tempB, aes(ymin = Lower, ymax = Upper, y = Mean), width = 0,
position = position_dodge(width = 0.6)) +
scale_shape_manual(values = c(19, 15)) +
scale_x_continuous(breaks = c(7:11), labels = c("7y", "8y", "9y", "10y", "adults")) +
labs(title = "Factor scores by age",
# subtitle = "Children (Study 2)\n",
x = "Age",
y = "Factor score") # 1000 by 500